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Abstract 

The influence of flexoelectric coupling on the internal structure of neutral domain walls in tetrag- 
onal phase of perovskite ferroelectrics is studied. The effect is shown to lower the symmetry of 
180-degree walls which are oblique with respect to the cubic crystallographic axes, while {100} and 
{110} walls stay "untouched". Being of the Ising type in the absence of the flexoelectric interac- 
tion, the oblique domain walls acquire a new polarization component with a structure qualitatively 
different from the classical Bloch-wall structure. In contrast to the Bloch-type walls, where the po- 
larization vector draws a helix on passing from one domain to the other, in the flexoeffect-affected 
wall, the polarization rotates in opposite directions on the two sides of the wall and passes through 
zero in its center. Since the resulting polarization profile is invariant upon inversion with respect 
to the wall center it does not brake the wall symmetry in contrast to the classical Bloch-type 
walls. The flexoelectric coupling lower the domain wall energy and gives rise to its additional 
anisotropy that is comparable to that conditioned by the elastic anisotropy. The atomic order- 
of-magnitude estimates shows that the new polarization component Pi may be comparable with 
spontaneous polarization P s , thus suggesting that, in general, the flexoelectric coupling should 
be mandatory included in domain wall simulations in ferroelectrics. Calculations performed for 
barium titanate yields the maximal value of the Pi , which is much smaller than that of the spon- 
taneous polarization. This smallness is attributed to an anomalously small value of a component 
of the "strain-polarization" elecrostictive tensor in this material. 
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I. INTRODUCTION 



In the light of miniaturization of electronic devices and achievements of nanotechnology, 
the question of functionality of domain walls in ferroelectrics is an exciting issue. The current 
limit for nano-scale engineering is of the order of 10s of nanometers. Having sizes typically 
of few nanometers, domain walls offer unique properties that are not exhibited in the bulk 
of a ferroic sample. For example, there is experimental evidence that twin domain walls in 
nonferroelectric CaTiC>3 possess spontaneous polarization 1 . Also ferroelectric properties are 
predicted in antiphase boundaries of otherwise nonferroelectric SrTiO^. As in the trend of 
miniaturization the next logical step is utilization of single domain wall as a functional ele- 
ment, the fundamental research of the domain wall's internal structure is highly demanded. 
If we consider for example 180-degree domain wall, which is the junction between two op- 
positely poled domains, its simplest profile contains only one polarization component (Ising 
wall). However it is possible that domain walls with additional polarization components are 
energetically favorable. A classic example is the Bloch wall, where an additional in-wall- 
plane polarization component arises, resulting in a helical polarization profile 3 . Perovskite 
crystals represent a class of materials with reach symmetry properties allowing a variety of 
wall structures. As mechanical and electric properties are strongly coupled in ferroelectrics, 
the domain wall structure is determined by both electric and elastic properties of a ma- 
terial, and taking into account the latter may radically affect the wall structure^. Up to 
now in the context of neutral ferroelectric domain walls the description of electromechanical 
coupling was restricted to the electrostrictive interaction. The electrostriction considerably 
influences the stability of Ising walls and introduces anisotropy of the wall energy^, but this 
effect does not introduce new features in the wall structure. Recent studies^ reveal a con- 
siderable impact of the generalized flexoelectricity (bilinear coupling between the strain and 
the order parameter gradient) on the wall structure in ferroics. In this paper we examine the 
effect of the flexoelectricity on electrically neutral ferroelectric domain walls as a function 
of wall orientation. We consider perovskite-type ferroelectrics in the tetragonal phase and 
perform numerical calculations for BaTiC>3 (BTO). We show that the flexoelectric effect has 
no impact on 90-degree walls and on 180-degree walls of the {100} and {110} orientation. 
At the same time for oblique 180-degree walls the effect leads to a wall with a new type of 
structure. 
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II. GINSBURG-LANDAU-DEVONSHIRE THEORY 

We consider perovskite material with m3m symmetry of the parent phase. The electric 
displacement field is defined as D = e^Ei + P, where £& is the background dielectric permit- 
tivity, E is the vector of electric field, P is the ferroelectric part of the polarization vector 
(hereafter we use the term of polarization as a shorthand). The D- field satisfies the Poisson 
equation: 

div~D = (1) 

The Gibbs free energy density expanded to sixth powers of polarization including gradient 
and flexoelectric terms is written as follows 5 . 



1 dP- dP k 

G = AijPiPj + BijklPiPjPkPl + CijklmnPiPjPkPlPmPn + 2^^ kl ~dl^~dx 



1 1 „ / dPu „ da. 



-QijklVijPkPl — TiSijkl&ijO'kl + nFijkl I Cij-l P) 



where Ajj = ai<%, = ^ < <%4z > +(ou - a 12 )g\f H , and Cy fcJm „ = sm < 5ij5 k i5 mn > 

+(0112 - ^) < MkL» > +( a m " fl H2 + *T)0i;Lm are the 2nd ' 4th ' and 6th order 
dielectric stiffness tensors. With <> we denote symmetrization with respect to interchange 

of indices: e.g. < 8ij6 k i >= \{8ij8 H + 8 ik 5ji + 8 u 5 jk ), 



D ijk i = Di 2 8ij8 k i + D e6 (8 ik 8ji + 8u8j k ) + (D n - D 12 - 2D m )g ( £ > kl (3) 

is the tensor controlling the correlation effects, are the components of mechanical stress. 
Hereafter the summation over repeating indices is implied, 8ij is the invariant Kronecker 
tensor, g^ kl and gf^ kimn are invariant tensors for the cubic symmetry. In the cubic crystallo- 
graphic axes the tensors g^ kl and g^ klmn have following structures: gjj kl — lifi — j — k — I 
and 9ijl = otherwise; gf) Umn = 1 if i=j = k = l = m = n and gf klmn = otherwise. 

Electrostriction tensor Qij k i, compliance tensor Sij k i and flexoelectric tensor Fij k i have 
same structures as D^f-f 

Qijki = Qvi8ij8 k i + ^-(8 ik 8ji + 8u8j k ) + (Q u - Q 12 - ^-)g\f k i (4) 
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Sijkl — S12&, 



'ijfal + —r(bikbji + SuSjk) + (su - S12 




(5) 



Fijki — F 12 Sij6ki H — — (SikSji + duSjk) + (Fu - 
From the Gibbs potential (j5J) one obtains equations of state: 



F12 — F§o)g, 



,(4) 
'ijkl 



(6) 




(7) 



For mechanical stresses we apply conditions of mechanical equilibrium: 



3xa 



(z,j = l-3) 



(8) 



In view of presence of stress gradient in the expression strain is defined via the variational 
derivation: 



III. STATEMENT OF THE PROBLEM FOR NEUTRAL DOMAIN WALLS 

We consider the material in the tetragonal phase, where the spatially homogeneous 
solution to the set of equations (C]),flE|) f° r mechanically free sample yields six equiva- 
lent domain states {P s ,0,0}, {— P s ,0,0}, {0,P S ,0}, etc. with spontaneous polarization 



the clapping angle) between the polarization vectors; between oppositely poled domains, 
180-degree walls are formed. Below we show that the flexoelectricity does not affect proper- 
ties of {100} and {110} walls. Condition of mechanical compatibility allows only one type 
of orientation for 90° walls, namely {110}. Thus 90° walls are not affected by flexoelectric 
coupling. That is why we consider only the 180-degree walls. 

A. 180-degree walls 




(9) 




. 90-degree walls separate domains with the angle of 90° (to within 



Electrically neutral 180-degree walls are parallel to the P s -vector. We characterize 
the orientation of the wall by the angle a between the wall normal and the OX 3 c cubic 
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FIG. 1. Orientation of the 180-degree domain wall and reference frames used. 

crystallographic direction as shown in Fig. [TJ Calculations are performed in the reference 
frame (0X l5 0X 2 , 0X 3 ) shown in Fig. [H which is associated with the wall. We consider a 
one-dimensional (ID) problem with the polarization vector P and mechanical stresses tensor 
<Jij being dependent only on the coordinate £3 normal to the plane of the wall. We neglect 
the polarization component normal to the wall as suppressed by the strong depolarizing 
field, so that only Pi and P2 components are allowed. In the new reference frame the Gibbs 
energy (J2]) reads 



G = ai {Pi + Pi) + a n P? + a 22 (a) P 2 4 + a 12 P?P2 + a in P? + a 222 (a) P 2 6 + (10) 
+a 112 P 4 P 2 2 + a 122 (a)P 2 4 P x + 



Dm 




2 + D u (a) 




2 


ydx 3 j 







-Qii^iPi 2 - Q22 (a) a 2 Pl - Q m a & P l P 2 - Q 12 (t^Pf + <r 2 P*) - 
~2 ( s n a i + s 22 (") °2 ) - s 12 a!a 2 - -s 66 cr 2 + -P 24 (a) (°"2^ - P 2^)- 

Here we omitted the terms that are null at 03 = 04 = a 5 = 0, since as we show below in 
Subsect. IIIIBl these stress components do not appear in the one- dimensional case. The 
designations used are: 

( \ 2<3n _ 0,12 ■ 2 m s 
0-22 («) = an sin (2a) , (11) 
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Q22(n) = On -sin 2 (2a) ('% _ QnQll) , (12) 



s 22 (a) = s u + sin 2 (2a) l^j- - Sn - Sl2 ) . ( 13) 

sin f 4d ) 

F 2223 (a) = F 24 (a) = (F 66 - F n + F 12 ) = F a sin (4a) (14) 



/ \ — ^112 9 / \ / \ 2(Zi^2 — ^123 9 / / \ 

«222 (a) = am sin (2a) , a i22 (a) = an 2 sin (2a) (15) 



D AA (a) = D 66 + sin 2 (2a) ( 11 g 12 - P 66 J (16) 

The following Voigt's (matrix) notations are used: 

an = ai, aim = an, 6an 22 = ai 2 , -Dim = -Dii, P1122 = P12, D\ 2 \ 2 = D^, Qim = Qu, 
Q1122 = Q12, 4Qi 2 i 2 = Sun = Sn, Sn 22 = Si 2 , 4si 2 i 2 = Sqq, Fim = F11, F1122 = F12, 

2^1212 = -^66 -~ 

From ffTUl) one obtains the equations of state for the polarization components depending 
only on x 3 in the form: 



2aiPi + 4a n P? + 2a 12 P^P 1 + 6a lu P? + Aa ll2 PfP 2 2 + 2a i22 (a) P X P 2 4 - (17a) 

<9 2 P 

-^e^ - 2 (q u <7i + g i2 a 2 ) Pi - g 66 a 6 p 2 = 

2a x P 2 + 4a 22 (a) P| + 2a 12 P 2 P 2 + 6a 222 (a) P| + 2a 112 P 1 4 P 2 + Aa 122 (a) P 2 P 2 3 - (17b) 

<9 2 P <9er 
-P44 (a) - 2 (Q 22 (a) a 2 + Q 12 a x ) P 2 - Q 66 <x 6 Pi ~ F 2i (a) ^ = 

The boundary conditions for the polarization far from the wall are 

Pi (x 3 ->• -00) = -P s , P 1 (x 3 ^oo) = P s , P 2 (x 3 ->• ±00) = 0. (18) 



B. Elimination of mechanical variables 

We consider the bulk of domains to be mechanically free: 

crijixs -»• ±00) = 0.i,j = 1,2,3. (19) 



This implies, using Eqs. ([2]) and fl9]), the boundary conditions for the strain components: 

£ 11 = QllPs'i £ 22 = £ 12 = 0- (20) 

For our ID problem, the condition of mechanical equilibrium (jS]) reads da 3 /dx 3 = 0, 
8(74/8x3 = 0, 8a 5 /8x 3 = 0. In view of (fT5j) . this condition requires that cr 3 = cr 4 = <7 5 = 
everywhere. The ID character of the problem, also enables us to rewrite, the Saint- Venant 
compatibility relationships 

ZikiZjmn {8 2 e ln /dx k dx m ) = 0, (21) 
[e-ijk is the Levi-Civita symbol) as: 

(feu/dxl = d 2 e 12 /dxl = d 2 e 22 /dxl = 0. (22) 
The solution to Eq. (|22|) with boundary conditions (fT9l) is: 

eifo) = Qii^l; e 2 {x 3 ) = Q 12 P$; e 6 {x 3 ) = (23) 

Note that it is the only possible one- dimensional solution for the elastic problem. The 
applicability of this solution to a stress-free finite sample is equivalent to the applicability of 
a one- dimensional model to a parallel plate capacitor. By applying this we neglect the fringe 
elastic fields at the contact of the domain wall with the surface, which is permissible when 
the dimensions of the sample are much larger than the thickness of the domain wall. Note 
that same "partially clamped" elastic conditions are usually applied for the description of 
mechanical stresses in a thin ferroelectric film on a substrate^. Under this ansatz we solve 
the system of equations and obtain expressions for the nonzero elastic stress components 
in the form: 



(Txi = 0i 



-F 24 (a) s 12 (dP 2 /dx 3 ) + (Pj - Pf) (Q u s 22 (a) - Q 12 s l2 ) 
+P 2 {Q22 (a) S12 - Qi 2 s 22 (a)) 

s 22 (a) Sn - s\ 2 

F24 (a) s n (0P 2 /dx 3 ) + (Q 12 s n - Q n s 12 ) (Pi - P?) 
+P| (-S11Q22 (a) + Q12S12) 



0"22 = 0"2 



s 22 (a) s n - s\ 2 



<7i2 = 0"6 



Q 



Mi 



S66 



Pi (x 3 ) P 2 (x 3 ) 



(24a) 

(24b) 
(24c) 



Eqs. ([21]) and (ITT)) form a full set of equations to define the polarization profile. 



IV. ANALYSIS OF EQUATIONS AND ESTIMATES 

First, as a benchmark, let us consider the set (JT7J) and ( I24j) without the flexoelectric 
effect (with F 24 (a) set to zero). One can check that at P 2 4(a) = the set (ITTj) and has 
single-component solution, as Eq. (If 7b j) can be satisfied with P 2 = 0. This Ising solution 



may be either stable or unstable^; in the latter case Bloch wall profile is observed with P 2 
being an even function of x%. We consider the first case where the Ising profile is stable if 
the flexoelectric effect is neglected. 

Since F 2 ^(a) is proportional to sin(4a) (ITT1) . the flexoelectric effect does not reveal itself 
for {100} (a = 0) and {110} (a = 7r/4) wall orientations. Hence, these walls remain Ising 
with the flexoelectric effect "switched on". For all the other wall orientations (hereafter 
termed as oblique walls), the flexoelectric coupling inevitably leads to the appearance of an 
additional polarization component. Indeed, substitution of ( )24|) into Eq. (117b[) introduces 
an additional coupling between the polarization components. As a result the one-component 
solution for the polarization profile (with P 2 = 0) is no longer available. 

It is instructive to note that {110} - oriented 90° domain walls do not "feel" the flex- 
oelectric coupling by the same reason. The statement of the problem for 90° walls is the 
same, but with boundary conditions 



Pi (x 3 ->• ±oo) = 0, P 2 (x 3 -> -oo) = -P./V2, P 2 (x 3 ^oo) = P s /V2 (25) 

instead of (JTHJ and with P 3 (x 3 ) = P 8 /y/2. One can carry out the same analysis as we 
have done for 180°-walls and check that flexoelectricity produce no nonzero terms in the 
equations of state for the polarization vector. Thus the flexoelectric effect does not affect 
either 90°-walls or {100} and {110}-oriented 180° walls. In contrast, oblique 180° walls 
inevitably acquire the new polarization component P 2 . 



A. Estimation of magnitude of fiexoelectric-effect-induced P 2 component 



To roughly estimate the magnitude of the P 2 component induced by the flexoelectric 



effect, we consider the linearized Eq. (117bl) . where we use for P 1 the profile of Ising wall, 



not perturbed by flexoelectric coupling. This approximation is fully valid at the vicinity of 
{100} and {110} wall orientations, where flexoelectric coupling is small. Neglecting P2 (273) 
with respect to Pi (x 3 ) in Eq.flT7aJ we obtain the standard wall profile^ 



Pi = Pi sinh 2 (x 3 /Rc) (cosh 2 (x 3 /R c ) + A) 1 (26a) 
Rc = V P§ (2a u + 6 6a m P|) (26b) 
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2oi 1 1 Pa 

A = — ? 26c 

011 = " n + 2s 22 (a)s u -4 2 (26d) 

Linearizing the set of Eqs. ( 124|) and ( 11 7b 11 with respect to P2 we obtain: 



(25! + 2a 12 P 2 + 2a 112 P 1 4 )P 2 - D u (a) ^ = / (a) ^ ^ (27a) 

Oi = ai + 6P 2 (27b) 

2i2 = Oi2 + 77~~~~ - © (27c) 

^66 

= g i2 (gig ~ Q11Q22 (a)) - Q12 (gaggii ~ Q11S22 (a)) ^ 27d ^ 



2s 22 (a) sn - s 



P44 (a) = P44 («) + F24 | Q) SU 2 (27e) 
s 22 (a) sn - s( 2 

f, \ (<5l2SH — QllSl2) „ , \ / Q „PN 

/ « = r~\ — F 24 a 27f 

S22 (a) Sn - s{ 2 

with Pi coming from ( 1261) . 

Pi being an odd function, the symmetry of Eq. (12T|) allows an odd solution for P2(xs) 
component. As we confirm by the numerical calculations below, the odd solution is stable, 
meaning in particular P 2 = 0; dP^jdxj, 7^ at £3 = 0. Hence, in the vicinity of £3 = 0, 
the main contribution to the Landau energy is due to the gradient term, which allows us to 
derive an approximate solution by neglecting the term linear with respect to P2. We also 
set A — > in Eq. (l30l for the simplicity. On simplifying thus, in the vicinity of £3 = Eq. 
(j27p transforms into 

»«'">S^ ,, " ) t (28) 
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The first integral of Eq. ( 128]) is 



P 44 (a) = / (a) ( (A 2 - Pj) + Co) • (29) 

We set Cq = to prevent the linear increase of P2 at 23 — > ±00. Taking into account that 
-^2(0) = 0, we obtain the solution in the form: 

P 2 (23) = -f^- T (if - PJ) dx « -^M-P|P C tanh f f ) (30) 

The solution (|30|) does not satisfy the boundary conditions P2 (23 — )■ ±00) = 0. It means 
that the other term that we do not take into account in fl28|) is responsible for the decay of 
P2 in the domains. However, the value P m = P2 (23 = 00) from (130]) may be used for the 
estimation of the amplitude of P 2 . Reverting to the initial designations we obtain P m in the 
form: 

Pm = 7 ^ ^ l \ ~ \ Ql2 ~ Qll-^-) P S R c (31) 

P44 (a) (s 22 (a) P 24 (a) 2 V Sl1 / 

We rewrite it in dimensionless form: 



£>44(a) 522(a) 

"11 



%£ = _^_2i^(^-£^ ) .r(a) (32a) 
P s VPeeSn V2aiiSn VQ11 s u 

4a)P> 66 /P> 44 (a) 

1-* + sin 2 (2a) (j*- + 



T (a) = -, sin(4c)D 66 /P 44 (a) ^ 




^) Wl-J + sin 2 (2a) + 



Let us analyze expression ( 1321) . The factor ^ D ^° sn is °f the order of unity according to 
atomic estimates. The applicability of atomic estimates for the evaluation of the electrostric- 
tive tensor in perovskite ferroelectrics is supported by experimental evidence of rather strong 
flexoelectric coupling in these materials 10 . The factor y|= , which is responsible for the 
sound velocity change near the ferroelectric phase transition is also of the order of unity in 
ferroelectrics with strong electromechanical coupling^. Landau theory does not require that 
the dimensionless factors — and T (a) are small compared to unity so that for rough 
estimation they may be taken as unity. Thus from formula ( 132]) the flexoelectric - driven 
component P2 is expected to be of the same order than P s . 
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It is instructive to specify the above estimates for BTO crystal, for which, further in the 
paper, we will present numerical simulations for the polarization profile of the domain bound- 
ary. For BTO the factor g^j — |^ appears to about 1/15. The smallness of this factor might 
be considered as purely accidental in view of seemingly different physical phenomena behind 
the elasticity and electrostriction. However, this smallness is of purely electromechanical na- 
ture. To see this we rewrite the factor QnjQw ~ Snjsw in terms of the "strain-polarization" 
electrostriction tensor linked with Qij by the relationship Qy = Sjiqu. The tensor can 
be considered as a primary material parameter since it is directly controlled by the lattice 
mechanics of the crystal. Being interested in perovskites where the Poisson ratio -s 12 /sn is 
typically about 0.3, we readily find Qn/Qn — s^/sii — 912/911 (0.5/(1 — Q.Q(qiij <iw))- Thus 
we see the smallness of the factor Qn/Qw ~ sw/sn is conditioned by the fact that, in metal- 
oxide ferroelectric perovskites, the ratio qyijqw is typically very small compared to unity. In 
particular, for SrTi0 3 q 12 /qw = -0.086^, in Pb(Zr0 3 )^ x (Ti0 3 ) x (PZT) \q 12 /q u \ < 0.05 
for x G (0.6, 1) (gn, qi% recalculated using Refs.— ,— ) . However there are materials where 
qi2 is not small with respect to qn. For example in KNbO^ qyijqw = —0.37, in tetragonal 
PZT near the morphotropic boudary qyijqw ~ —0.5 ( 13 , 14 ). There is no data available 
on the flexoelectric tensor coefficients in the latter materials, but one can expect there 
relatively large domain-wall energy anisotropy and flexoelectric-effect-induced polarization 
components due to the elevated qyijqw ratio. 

For BTO the estimate for the amplitude of second polarization components is also affected 
by the exceptionally high anisotropy of the correlation energy {D\\jD m ~ 25). This occurs 
via the factor D ^ a ^ in the estimate ( 132]) . Thus we conclude that as being additionally 
affected by two small factors, the maximal value of the second polarization component in 
180° walls in BTO is expected to be two orders of magnitude smaller than P s . We would 
like to stress that there is no reason to expect anomalously small values of this component 
in ferroelectrics in general. 



V. NUMERICAL RESULTS AS APPLIED TO BARIUM TITANATE 

To be more specific, we analyze the impact of the flexoelectric coupling on the structure 
of a domain wall in ferroelectrics for the case of the tetragonal phase in classical perovskite 
ferroelectric BaTiOa at room temperature. Since the problem is not analytically tractable, 
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we do it numerically using the thermodynamic parameters of this crystal. Since the experi- 
mental values of components of the flexocoupling tensor F^m are not currently available, we 
use in our calculations the F^i tensor evaluated from the data of ab initio calculations for 
Bao.sSro.sTiOs (BST) from Ref.— . In case where the results are found critical to the exact 
value of this tensor, we vary these to cover possible situations. The parameters used for the 
simulations are listed in the Table 1. 

The profiles for Pi and P2 polarization components obtained from a numerical solution 
to Eqs. ( ITT)) and (!!MI) are shown in Fig. [2j We use (fT8l) as the boundary conditions with an 
additional condition of vanishing of spatial derivatives of all variables at the infinity. The Pi 
profile is perfectly described by formula ( )26|) . the difference conditioned by the flexoelectric 
effect is within the line width of the plot for any angle a. The maximal value of P2 is as 
small as 6 • 10~ 2 fiC / cm 2 which is ~ Ps/300. One can see from the Fig. [2] that the width 
of the domain wall with respect to P2 is few times larger than with respect to Pi. This is 
a consequence of the small ratio of correlation lengths 454 = / D66 which results from 
Du/D m ~ 25. In the limit P2 <C Ps and r(P{) <C r(P 2 ) an approximate analytical solution 
for the wall profile can be also developed, which is given in Appendix |A] Dependence of 
the maximal polarization P2 (x^) on the angle a is shown in Fig. |3j The anisotropy may 
be understood from expression ( |3T|) as an interplay between the angular dependence of the 
flexoelectric factor F 2 4 («) and that of the gradient-energy factor .D 44 (a). 

Inspection of the results of the simulations reveals that the flexoelectric coupling can 
lead to the formation of a polarization profile which has never been obtained to ferroelectric 
domain walls. The polarization rotates in the boundary, but the latter does not acquire chi- 
rality. The senses of the polarization rotation on the two sides of the boundary are opposite, 
the polarization passing though zero at its center. As the result the polarization profile of 
the modified wall remains invariant with respect to the inversion about the wall center, in 
contrast to the Bloch-type wall. Following the terminology introduced by Houchmandzadeh 
et al^ such wall can be classified as bichiral. One should note that the rough approximate 
solution fl30|) qualitatively correctly reproduced the bichirality of the wall. 
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FIG. 2. Pi(xs) and ^2(^3) - profiles across 180° domain wall calculated for the parameters of 
BaTi03 at room temperature given in Table 1. The Pi(xs) -profile is practically independent 
of the angle a between the wall and a cubic crystallographic direction. The P2(x$) - profile is 
calculated for a = ir/8 



A. Free energy calculations 



In this subsection we study the anisotropy of DW energy induced by flexoelectric effect. 
The expression for the Free energy density $ may be obtained from Gibbs potential G by 
Legendre transformation $ = G + <JiSi. The wall energy per unit area Ew is then given 
by the integral Ew = /_^°°(^(^3) — $00)^3 where $oo is the energy density at 23 — > ±00. 
From (|T0|) with the aide of (|23|) and (|24"|) we derive the following expression for the wall free 
energy: 
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Angle a (rad) 



FIG. 3. Maximal absolute value of of the second polarization component in the wall as a function 
of the angle a between it and a cubic crystallographic direction calculated for the parameters of 
BaTiO-3 at room temperature given in Table 1. (a) - plot in one quadrant; (b) - polar plot. 



E w = 



poo 

I {a, (P 2 - P 2 + Pi) + a n (P x 4 - P 4 ) + a 22 (a) P? + a 12 P?Pi + (33) 

— oo 

+o m (Pf - Pi) + a 222 (a) P 2 6 + a 112 P 4 P 2 2 + a 122 (a) P 2 4 P 2 + 
2 \dx 3 J 2\ s 22 {a)s u -s( 2 J \dx 3 J 

12 s ll — Qll s 12 fT-,2 t-.2\ dP"2 , {Qll (Pf ~ Pf) ~ Ql2P 2 ) 



" \ ^ A / \ U Z2 "11 "12/ \ / 

LZ? /x Q\2S\\ - Q\\S\2 / p 2 p2\ <9p2 (Qll - A 2 ) - Ql2P2 2 ) 2 
+-^24 (a) -T~\ 2~ l-Ps _ M J a 1 o 

S22 (a) sn - sf 2 cte 3 2s n 

((Q12S11 - giifig) (Pj ~ Pi 2 ) + P 2 2 (-5nQ 22 (a) + Qi 2 si 2 )) 2 + g|^ p2p2} ^ 
2su (s 22 (a) su - s? 2 ) 2s 66 1 2 



Using the coordinate dependences obtained numerically and this relationship, we calculate 
the energy of the bichiral wall as a function of the angle a plotted in the Fig. |5j It is seen 
from this figure that the anisotropy of the wall energy is extremely weak. However, it is 
worth looking closer the fine stricture of the anisotropic part of the wall energy, specifically 
to decompose it to the parts independent of and dependent on the flexoelectric coupling E$ A 
and Ey^ A , respectively. Such decomposition is presented in FigElfbr the case where for the 
flexoelectric coefficients of BTO we use those obtained by ab initio calculations for BST. As 
we can see there is an essential difference between the E^) A and Eyp A anglular dependences. 
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FIG. 4. Schematic of the structures of neutral 180-degree domain walls addressed in the paper, (a) 
- the Ising type structure occurring when the wall is normal to the cubic crystallographic directions 
or/and when the flexoelectric coupling is isotropic or neglected; (b) - the bichiral structure occurring 
for the oblique orientation of the wall provided that the flexoelectric coupling is anisotropic. 

Ey$A IS minimal at a = 0, maximal at a = ir/A and monotonic in the interval (0, vr/4), while 
E)^ A is maximal in both a = and a = tt/4, and has a minimum in an intermediate point, 
around a = 7r/10. Remarkably, for small angles a the both contributions are comparable, 
the flexoelectric one being a bit smaller. However if we take the flexoelectric factor F a two 
times as large as that for BST, flexoelectric contribution becomes dominating at small a. 
This leads to qualitative change of the anisotropy of the wall energy: flexoelectric coupling 
results in splitting each of the energy minima (at a = and a = vr/2) into two, see Fig. UJ 
The results of the numerical calculations can be elucidated using some analytical rela- 
tionships. First, the contribution E^) A can be readily presented in a simple form 

EwA -^{O^ TTj s 22 (a) s n - s\ 2 L ^ ~ Pl > dX3 (34) 
which is consistent with the results by Dvorak and Janoveo^. As was recognized by these 
authors (and clear from this expression), the angular dependence of this contribution is 
conditioned by that of the elastic compliance. The flexoelectricity-conditioned contribution 
can be evaluated taking into account that in the case of BTO P 2 P s and r(Pi) -C r(P 2 ). It 
enables us to keep among P 2 -containing terms in f )33|) only dP 2 /dx s (one can check that terms 
containing P| are smaller by factor jy^^) while using (|29|) as an approximate relationship 
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FIG. 5. Energy of bichiral wall as a function of the angle a between it and a cubic crystallographic 
direction. Parameters used in the calculations are listed in Table 1 at room temperature T=293 K. 

for dPijdxz to get: 



<F) „ 2F 24 («) s n 



E WA « -E^at; r—7 7 s 7T- ( 35 ) 

D44 (a) (s 22 (a) sii - sf 



It is clear from this relationship that the smallness of both anisotropic contributions to the 
wall energy is controlled by the square of the same factor Q12/Q11 — Si2/ S n — 1/15, which 
was already recognized responsible for the smallness of the second polarization component. 
As for the shape of the angular dependence of the wall energy, it is conditioned by an 
interplay among the elastic anisotropy, anisotropy of the correlation energy, and that of the 
flexoelectric coupling. 



VI. CONCLUSIONS 

The presented analysis demonstrates that the flexoelectric coupling once it is anisotropic 
(at F a = Fu — F\i — F 66 = all effects addressed in the paper disappear) results in the 
formation of 180-degree ferroelectric domain walls where the polarization rotates but the 
wall does not acquire chirality. Following the terminology introduced by Houchmandzadeh 
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FIG. 6. Comparison of the anisotropic contributions to the energy of bichiral wall as functions of the 

(s) 

angle a between it and a cubic crystallographic direction. EyJ A is the contribution solely controlled 
by the elastic anisotropy. E^y A is the contribution due to the anisotropy of the flexoelectric 
coupling. The calculations are done for the parameters listed in Table 1 at room temperature 
T=293 K. 

et al 1 ™ such wall can be classified as bichiral. In contact to the Bloch type walls, the 
appearance of the second polarization component in bichiral walls dose not brake the wall 
symmetry with respect to the spatial inversion. In addition, depending on the value of F a , the 
flexoelectric coupling can lead to the doubling of the number of the energetically favorable 
orientations of the walls. Order-of-magnitude estimates show that, in general, the effects 
driven by the flexoelectric coupling can be appreciable. However, for the thermodynamic 
parameters of barium titanate, the calculations performed show that the amplitude of the 
second component is expected to be smaller than one percent of the spontaneous polarization 
while the modulation of the wall energy (as a function of its orientation) is found to be yet 
smaller. The numerical smallness of these effects is shown to be mainly conditioned by that 
of the qi2 component of the electrostriction " strain-polarization" tensor, which is typical for 
metal-oxide ferroelectric perovskites. We would like to stress that there is no reason to expect 
the flexoelectricity-induced features of 180-degree ferroelectric domain walls addressed in this 
paper to be anomalously small in ferroelectrics in general. 

18 




-0.3 0.3 

Angle a (rad) 



FIG. 7. Energy of bichiral domain wall as a function of its orientation for different values of 
the flexoelectric coupling calculated for the parameters of BTO (Table 1) at room temperature 
T=293. The curves labeled with 0, 1, and 2 correspond to the values of the anisotropic part of 
the flexoelectric coupling F a equal to 0, Fa BST \ and 2Fa BST ^ , where Fa BST ^ is the value of the 
flexoelectric coefficient for BST recalculated using the results of ab initio from Ref.— (see Table 1). 
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Appendix A: Approximate analytical solution for the wall profile 

Let us derive an analytical solution for the polarization profile in the approximation 
P 2 <C Pg] r(Pi) <C r(P 2 ). We solve Eq. (l2Tj) taking Pf — > P| everywhere except the 
derivative which is the driving force for the P 2 appearance. This approximation is 
valid far from the wall center where P\ ~ P s , and in the wall center where the gradient term 
is dominating over the term distorted by Pf — > Pg. 

The solution to the second order differential equation with constant coefficients, coordi- 
nate dependent inhomogeneity and boundary conditions ( JTB1 could be written using Green 
function method 19 . Rewriting the equation (127]) for P 2 in the from 
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p_ r 2&£2 = ^24 (a) (Ql2Sll-QllS 12 ) d {Pi -Pi) 

2 T dx\ 2a 1 + {2a 12 + {Q 2 66 /s m ))Pl + 2a n2 Pl s 22 {a) s n - s\ 2 dx 3 1 } 



where 



-D44 («) 



2 ai + (2a 12 + (Q 2 Js m )) P| + 2a 112 P| 



one could readily 6 find that 



R 2 c -4r 2 \ "V P, 



exp ( — — ] — exp ( — ■ 



(A2) 



ftM 'f^f-^ M^ . (A3) 



P24 (oQ (Qi2gn - giifig) 

2ai + (2o 12 + (Qge/see)) P| + 2a 112 Pj 

s 22 (a) Sn s\ 2 

Using an approximation for Pi profile (1261) 

Pi « P s (1 - exp (- |x 3 | /Rc)) sign(x 3 ) (A4) 
we obtained from Eq.(15b) the following expression: 

P( T \~ -f(a)P 2 p(x 3 ) 

^ W ~ 2ai + (2a 12 + {Q 2 Js m )) P| + 2a 112 Pj 1 ° &j 

pM = (exp (-gj - exp (-^) J + (A5b) 

2P C / / 2x 3 \ ( x 3 



Comparison of the approximate analytical solution (1A5|) with numerical calculations 
based on the coupled equations (ITT)) is shown in Fig. |HJ The difference between the ap- 
proximate analytical expression flA5j) (dashed curve) and numerical simulation (solid curve) 
is of the order of several percent. Thus, for P 2 (x 3 ) analytical expression (1A5j) works with 
sufficient accuracy. 
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TABLE I. Free energy coefficients for bulk ferroelectric BaTiC>3 (from Refs.—,— ,—,—). 



parameters 


values 


Refs. and Notes 


ai (C- 2 mJ) 


a x = 3.34(T - 381) x 10 5 


Ref^ 


a { j (C" 4 m 5 J) 


an ~ 4.69(T - 393) x 10 6 - 2.02 x 10 8 
012 ~ 3.230 x 10 s 


Ref.i^ 


<Hjk (C" 6 m 9 J) 


am ~ -5.52(T - 393) x 10 7 + 2.76 x 10 9 
a 112 = 4.47 x 10 9 
ai2 3 = 4.91 x 10 9 






D$ = 5.1 x 10- 10 
D$ ~ -0.2 x 10- 10 
~ 0.2 x 10- 10 


Ref.-i 7 - 


Q t3 {C- 2 m A ) 


Qu = 0.11, Q12 ~ -0.043, Q 66 = 0.059 


Rcf. 7 


8ij{lQr vl Pcr x ) 


sn = 8.3,si2 ^ -2.7, s 66 = 9.24 


Ref.I 




/u = 5.12, / 12 = 3.32, / 66 = 0.045 


Ref.i^ 




F n = 2.46, F12 = 0.48, F m = 0.05 

P{ BST ) _ _ ^ _j_ _1_93 


Recalculated us- 
ing relationship 

Forf = f/3-yS/S-y', f a .y 

taken from RefJ^ 


(C- 2 m 3 J) 


= D u = 3.52 x 10~ 10 
= -°12 = -1-24 x lO" 10 
D$ = D m = 0.2 x 10- 10 


recalculated using 

n (°") nM f zp 
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